import ipywidgets
%matplotlib inline
import matplotlib.pyplot as plt
import pandas as pd
import numpy as np
from tqdm.auto import tqdm
import igraph as ig
import xnet as xn
import importlib
import math
%config InlineBackend.figure_format = 'retina'
from collections import Counter
importlib.reload(xn)
dataSuffix = "2019" #2019 = "", 2020 = "2020"
maxClassSize = 69;
applyFilterMasks = False;
sizeProperty = "CLASS_SIZE_ENROLLED"
# sizeProperty = "CLASS_CAPACITY"
#Loading networks
gstudents = xn.xnet2igraph("Networks/student2student_classes_intersection%s.xnet"%dataSuffix);
# xn.xnet2igraph("Networks/student2student_classes_intersection_withHousing.xnet");
gclasses = xn.xnet2igraph("Networks/classes2classes_students_intersection%s.xnet"%dataSuffix);
#Loading *2* data
students2Classes = [];
classes2Students = [];
# with open("Data/students2Classes.txt","r") as fd:
# for line in fd:
# line = line.strip();
# if(line):
# entries = line.split("\t");
# students2Classes.append([int(entry) for entry in entries]);
# else:
# students2Classes.append([]);
with open("Data/classes2Students%s.txt"%dataSuffix,"r") as fd:
for line in fd:
line = line.strip();
if(line):
entries = line.split("\t");
classes2Students.append([int(entry) for entry in entries]);
else:
classes2Students.append([]);
gclasses.vs["Students"] = classes2Students;
#Filter classes by size and extra metadata:
maskSize = (np.array(gclasses.vs[sizeProperty])>maxClassSize)
mask = maskSize
if(applyFilterMasks):
maskMulti = [entry=="YLEC" for entry in gclasses.vs["MULTI_COMPONENTS_FLAG"]]
maskCLab = [entry.startswith("Y") for entry in gclasses.vs["COMPUTER_LAB_IND"]]
mask= mask+maskMulti+maskCLab;
vertices_toDelete = np.where(mask)[0]
gclasses.delete_vertices(vertices_toDelete);
gclasses.vcount()
#Calculate some properties
gclasses.vs["Strength"] = gclasses.strength(weights="weight")
gclasses.vs["Degree"] = gclasses.degree()
display("calculating closeness...");
gclasses.es["distance"] = 1.0/np.array(gclasses.es["weight"])
# gclasses.vs["Closeness"] = gclasses.closeness(weights="distance")
display("calculating betweenness...");
gclasses.vs["Betweenness"] = gclasses.betweenness(weights="distance",directed=False)
#Random for baseline
gclasses.vs["Random"] = np.random.random(gclasses.vcount())
gclasses.vs["kcore"]= gclasses.coreness()
classes2Students = gclasses.vs["Students"]
%load_ext Cython
%%cython
#--annotate
cimport numpy as np
import cython
cimport cython
ctypedef np.float64_t DTYPE_t
from libc.stdlib cimport rand, RAND_MAX
from libc.math cimport exp
cdef extern from "stdlib.h":
double drand48()
void srand48(long int seedval)
@cython.boundscheck(False)
@cython.wraparound(False)
@cython.nonecheck(False)
@cython.cdivision(True)
cdef double randomUniform():
cdef double r = rand()
return r / RAND_MAX
def setRandomSeed(long seed):
srand48(seed);
@cython.boundscheck(False)
@cython.wraparound(False)
@cython.nonecheck(False)
@cython.cdivision(True)
def updateClass(long[::1] studentStates, double[::1] lastUpdate,long[::1] classesGroupStudents, double gamma,double beta,double beginTime, double endTime):
# localStates = studentStates[classesGroupStudents];
# nS = localStates==-1
cdef long nI = 0;
cdef long studentIndex;
cdef double durationSinceUpdate;
cdef double flip;
cdef double rRate;
cdef double duration;
cdef long oldState;
for inGroupIndex in range(classesGroupStudents.shape[0]):
studentIndex = classesGroupStudents[inGroupIndex];
if(studentStates[studentIndex] == 0):
nI+=1;
durationSinceUpdate = endTime - lastUpdate[studentIndex]
if (durationSinceUpdate>0):
flip = drand48();
rRate = 1.0 - exp(-gamma*durationSinceUpdate)
lastUpdate[studentIndex] = beginTime;
if flip < rRate:
studentStates[studentIndex] = 1;
nI -= 1
for inGroupIndex in range(classesGroupStudents.shape[0]):
studentIndex = classesGroupStudents[inGroupIndex];
duration = endTime - beginTime
oldState = studentStates[studentIndex]
if (duration>0):
if (oldState == -1):
flip = drand48()
rRate = 1.0 - exp(-beta*duration*nI)
if flip < rRate:
lastUpdate[studentIndex] = endTime
studentStates[studentIndex] = 0
# nI += 1 #enable for sync
# if (oldState == -1):
# flip = randomUniform()
# rRate = 1.0 - exp(-gamma*duration)
# lastUpdate = endTime
# if flip < rRate:
# studentStates[studentIndex] = 1
# nI -= 1
@cython.boundscheck(False)
@cython.wraparound(False)
@cython.nonecheck(False)
@cython.cdivision(True)
def updateRecovery(long[::1] studentStates, double[::1] lastUpdate, double gamma, double currentTime):
# localStates = studentStates[classesGroupStudents];
# nS = localStates==-1
cdef long studentIndex;
cdef double durationSinceUpdate;
cdef double flip;
cdef double rRate;
for studentIndex in range(studentStates.shape[0]):
if(studentStates[studentIndex] == 0):
durationSinceUpdate = currentTime - lastUpdate[studentIndex]
if (durationSinceUpdate>0):
flip = drand48();
rRate = 1.0 - exp(-gamma*durationSinceUpdate)
lastUpdate[studentIndex] = currentTime;
if flip < rRate:
studentStates[studentIndex] = 1;
@cython.boundscheck(False)
@cython.wraparound(False)
@cython.nonecheck(False)
@cython.cdivision(True)
def countSIR(long[::1] studentStates):
# localStates = studentStates[classesGroupStudents];
# nS = localStates==-1
cdef long nS = 0;
cdef long nI = 0;
cdef long nR = 0;
for studentIndex in range(studentStates.shape[0]):
if(studentStates[studentIndex] == -1):
nS+=1;
elif(studentStates[studentIndex] == 0):
nI+=1;
elif(studentStates[studentIndex] == 1):
nR+=1;
return (nS,nI,nR);
def generateClassesSchedule(p=0.1,propertyName = "Betweenness Centrality"):
startTimes = gclasses.vs["CLS_START_TIME"]
endTimes = gclasses.vs["CLS_END_TIME"];
meetDays = gclasses.vs["CLS_MEETING_DAYS"];
# classSizes = gclasses.vs["Size"];
rankingProperty = np.array(gclasses.vs[propertyName]);
classesByDay = {letter:[] for letter in "MTWRFSN"};
validClasses = np.where(np.array([d!="nan" for d in meetDays]))[0]
rankingIndices = rankingProperty[validClasses].argsort()
maxRank = round((1-p)*len(validClasses))
selectedIndices = validClasses[rankingIndices][0:maxRank]
mask = np.ones(gclasses.vcount(), np.bool)
mask[selectedIndices] = 0
notSelectedIndices = np.arange(gclasses.vcount())[mask]
gclassesCopy = gclasses.copy();
gclassesCopy.delete_vertices(notSelectedIndices);
ggiantClusters = gclassesCopy.clusters()
numberComponents = len(ggiantClusters)
sizeLargestComponent = ggiantClusters.giant().vcount()
for classIndex in selectedIndices:
if(meetDays[classIndex]!="nan"):
for weekday in meetDays[classIndex]:
if(weekday=="D"):
for wlist in classesByDay:
classesByDay[wlist].append(classIndex);
else:
classesByDay[weekday].append(classIndex);
#reordering
studentsInClassesByDay = {};
startTimesOfClassesByDay = {};
endTimesOfClassesByDay = {};
for weekday in classesByDay:
classesByDay[weekday] = sorted(classesByDay[weekday], key=lambda index: startTimes[index])
studentsInClassesByDay[weekday] = [np.array(classes2Students[classesIndex],dtype=np.int) for classesIndex in classesByDay[weekday]];
# startTimesOfClassesByDay[weekday] = ["%s: %f"%(startTimes[i][-8:-6]+startTimes[i][-5:-3],float(startTimes[i][-8:-6])+float(startTimes[i][-5:-3])/60) for i in classesByDay[weekday]]
startTimesOfClassesByDay[weekday] = [float(startTimes[i][-8:-6])+float(startTimes[i][-5:-3])/60 for i in classesByDay[weekday]]
endTimesOfClassesByDay[weekday] = [float(endTimes[i][-8:-6])+float(endTimes[i][-5:-3])/60 for i in classesByDay[weekday]]
return classesByDay,studentsInClassesByDay,startTimesOfClassesByDay,endTimesOfClassesByDay,numberComponents,sizeLargestComponent
weekorder = "MTWRFSN"
def iterateModel(seed=100,gamma = 1.0/6.0/24.0/2.0/(5/7), beta = 0.33/24.0, startProbability = 0.00025, weeksCount = 12, trackBeforeClasses = False):
setRandomSeed(seed);
currentDay = 0;
studentStates = -np.ones(studentsCount,dtype=np.int)
lastUpdate = np.zeros(studentsCount)
studentStates[np.random.random(len(studentStates))<startProbability] = 0
studentHaveClasses = np.zeros(studentsCount,dtype=np.bool)
for week in weekorder:
for classesGroupIndex,classesGroupStudents in enumerate(studentsInClassesByDay[week]):
studentHaveClasses[classesGroupStudents]=True;
studentStates[~studentHaveClasses]=-2 #Invalid These students do not participate in any classes
times = [];
valuesOverTime = [];
if(not trackBeforeClasses): #First values
times.append(0.0);
valuesOverTime.append(countSIR(studentStates));
for weekIndex in range(weeksCount):
for week in weekorder:
# lastClassesStartTime=-10;
for classesGroupIndex,classesGroupStudents in enumerate(studentsInClassesByDay[week]):
startTime = currentDay*24+startTimesOfClassesByDay[week][classesGroupIndex];
endTime = currentDay*24+endTimesOfClassesByDay[week][classesGroupIndex];
if(trackBeforeClasses):
times.append(startTime/24.0);
valuesOverTime.append(countSIR(studentStates));
updateClass(studentStates, lastUpdate,classesGroupStudents , gamma,beta,startTime, endTime);
currentDay+=1;
updateRecovery(studentStates, lastUpdate, gamma, currentDay*24);
times.append(currentDay);
valuesOverTime.append(countSIR(studentStates));
nS,nI,nR = [np.array(entry) for entry in zip(*valuesOverTime)]
return times,nS,nI,nR;
# from collections import Counter
# studentProperty = gstudents.vs["STUDENT_LEVEL"]
# studentKeys = dict(Counter(studentProperty).most_common(8)).keys()
# studentDistribs = {key:[] for key in studentKeys}
# studentDistribs["Other"] = []
# validStudentProperty = [studentProperty[i] for i,valid in enumerate(studentHaveClasses) if valid];
# currentDistrib = dict(Counter(validStudentProperty))
# remainingStudents = len(validStudentProperty)
# for key in studentKeys:
# if(key in currentDistrib):
# keyCount = currentDistrib[key]
# studentDistribs[key].append(keyCount);
# remainingStudents-=keyCount;
# studentDistribs["Other"].append(remainingStudents);
# studentDistribs
import multiprocessing
from multiprocessing import Pool
from functools import partial
import random
import sys
from collections import Counter
studentProperty = gstudents.vs["STUDENT_LEVEL"]
gamma = 1.0/6.0/24.0/3.0#/24.0/2.0/(5/7)
beta = 0.33/24.0
startProbability = 0.00025
removedClassesRatio = 0.3
propertyName = "Betweenness"
realizationsCount = 500;
studentsCount = gstudents.vcount();
studentKeys = dict(Counter(studentProperty).most_common(8)).keys()
studentDistribs = {key:[] for key in studentKeys}
studentDistribs["Other"] = []
classesByDay,studentsInClassesByDay,startTimesOfClassesByDay,endTimesOfClassesByDay,numberComponents,sizeLargestComponent = generateClassesSchedule(removedClassesRatio,propertyName)
studentHaveClasses = np.zeros(studentsCount,dtype=np.bool)
for week in weekorder:
for classesGroupIndex,classesGroupStudents in enumerate(studentsInClassesByDay[week]):
studentHaveClasses[classesGroupStudents]=True;
validStudentProperty = [studentProperty[i] for i,valid in enumerate(studentHaveClasses) if valid];
currentDistrib = dict(Counter(validStudentProperty))
remainingStudents = len(validStudentProperty)
for key in studentKeys:
if(key in currentDistrib):
keyCount = currentDistrib[key]
studentDistribs[key].append(keyCount);
remainingStudents-=keyCount;
studentDistribs["Other"].append(remainingStudents);
num_processors = multiprocessing.cpu_count()
p=Pool(processes = num_processors)
allResults = p.map(partial(iterateModel,beta=beta,gamma=gamma,startProbability = startProbability),[np.random.randint(-sys.maxsize-1,sys.maxsize+1) for i in range(realizationsCount)])
p.terminate()
# allResults = []
# for i in tqdm(range(3)):
# allResults.append(iterateModel(0,seed=np.random.random()));
times,avgS,avgI,avgR = np.average(allResults,axis=0);
_,stdS,stdI,stdR = np.std(allResults,axis=0);
numStudents = np.average(avgS+avgI+avgR);
np.array(allResults)[:,2,:]
peakMagnitudes = np.max(np.array(allResults)[:,2,:],axis=1)
peakPositions = np.argmax(np.array(allResults)[:,2,:],axis=1)
avgPeakMagnitude = np.average(peakMagnitudes)
stdPeakMagnitude = np.std(peakMagnitudes)
totalInfected = (np.array(allResults)[:,2,-1]+np.array(allResults)[:,3,-1])
avgTotalInfected = np.average(totalInfected)
stdTotalInfected = np.std(totalInfected)
avgPeakPosition = np.average(peakPositions)
stdPeakPosition = np.std(peakPositions)
print("Num Students: %d"%numStudents)
print("Num Students2: %d"%np.sum(studentHaveClasses))
print("Total Infected: %f ± %f"%(avgTotalInfected,stdTotalInfected))
print("Peak Magnitude: %f ± %f"%(avgPeakMagnitude,stdPeakMagnitude))
print("Peak Position: %f ± %f"%(avgPeakPosition,stdPeakPosition))
print("Number of comp.: %d"%numberComponents)
print("Size giant comp.: %d"%sizeLargestComponent)
print("Distrib: "+str(studentDistribs))
plt.close("all")
plt.figure()
topS = avgS+stdS
lowS =avgS-stdS
lowS[lowS<0] = 0
topI = avgI+stdI
lowI =avgI-stdI
lowI[lowI<0] = 0
topR = avgR+stdR
lowR =avgR-stdR
lowR[lowR<0] = 0
plt.plot(times,avgS,label="S",color="#848484");
plt.fill_between(times,lowS, topS,alpha=0.25,color="#848484");
plt.plot(times,avgI,label="I",color="#AA403D");
plt.fill_between(times,lowI, topI,alpha=0.25,color="#AA403D");
plt.plot(times,avgR,label="R",color="#7EA0CA");
plt.fill_between(times,lowR, topR,alpha=0.25,color="#7EA0CA");
plt.legend()
plt.xlabel("Days")
plt.show()
#Over class removal ration
import multiprocessing
from multiprocessing import Pool
from functools import partial
import random
import sys
from collections import Counter
classesRatios = np.linspace(0,1,20,endpoint=False)
gamma = 1.0/6.0/24.0 #/2.0/(5/7)
beta = 0.33/24.0
startProbability = 0.0005
realizationsCount = 100;
studentProperty = gstudents.vs["STUDENT_LEVEL"]
# possibleProperties = ["Strength","Degree","Closeness","Betweenness","Size","kcore","Random"]
possibleProperties = ["Strength","Degree","Betweenness","Size","Random"]
numStudentsArrays={}
avgPeakMagnitudeArrays={}
stdPeakMagnitudeArrays={}
avgTotalInfectedArrays={}
stdTotalInfectedArrays={}
avgPeakPositionArrays={}
stdPeakPositionArrays={}
numberComponentsArrays={}
sizeLargestComponentArrays={}
studentDistribsByProperty={};
pbar = tqdm(total=len(possibleProperties)*len(classesRatios));
for propertyName in tqdm(possibleProperties):
studentsCount = gstudents.vcount();
studentKeys = dict(Counter(studentProperty).most_common(8)).keys()
studentDistribs = {key:[] for key in studentKeys}
studentDistribs["Other"] = []
numStudentsArray = [];
avgPeakMagnitudeArray = [];
stdPeakMagnitudeArray = [];
avgTotalInfectedArray = [];
stdTotalInfectedArray = [];
avgPeakPositionArray = [];
stdPeakPositionArray = [];
numberComponentsArray = [];
sizeLargestComponentArray = [];
for removedClassesRatio in tqdm(classesRatios):
pbar.update(1)
classesByDay,studentsInClassesByDay,startTimesOfClassesByDay,endTimesOfClassesByDay,numberComponents,sizeLargestComponent = generateClassesSchedule(float(removedClassesRatio),propertyName)
studentHaveClasses = np.zeros(studentsCount,dtype=np.bool)
for week in weekorder:
for classesGroupIndex,classesGroupStudents in enumerate(studentsInClassesByDay[week]):
studentHaveClasses[classesGroupStudents]=True;
validStudentProperty = [studentProperty[i] for i,valid in enumerate(studentHaveClasses) if valid];
currentDistrib = dict(Counter(validStudentProperty))
remainingStudents = len(validStudentProperty)
for key in studentKeys:
if(key in currentDistrib):
keyCount = currentDistrib[key]
studentDistribs[key].append(keyCount);
remainingStudents-=keyCount;
else:
studentDistribs[key].append(0);
studentDistribs["Other"].append(remainingStudents);
num_processors = multiprocessing.cpu_count()
p=Pool(processes = num_processors)
allResults = p.map(partial(iterateModel,beta=beta,gamma=gamma,startProbability = startProbability),[np.random.randint(-sys.maxsize-1,sys.maxsize+1) for i in range(realizationsCount)])
p.terminate()
# allResults = []
# for i in tqdm(range(3)):
# allResults.append(iterateModel(0,seed=np.random.random()));
times,avgS,avgI,avgR = np.average(allResults,axis=0);
_,stdS,stdI,stdR = np.std(allResults,axis=0);
numStudents = np.average(avgS+avgI+avgR);
peakMagnitudes = np.max(np.array(allResults)[:,2,:],axis=1)
peakPositions = np.argmax(np.array(allResults)[:,2,:],axis=1)
avgPeakMagnitude = np.average(peakMagnitudes)
stdPeakMagnitude = np.std(peakMagnitudes)
totalInfected = (np.array(allResults)[:,2,-1]+np.array(allResults)[:,3,-1])
avgTotalInfected = np.average(totalInfected)
stdTotalInfected = np.std(totalInfected)
avgPeakPosition = np.average(peakPositions)
stdPeakPosition = np.std(peakPositions)
numStudentsArray.append(numStudents)
avgPeakMagnitudeArray.append(avgPeakMagnitude)
stdPeakMagnitudeArray.append(stdPeakMagnitude)
avgTotalInfectedArray.append(avgTotalInfected)
stdTotalInfectedArray.append(stdTotalInfected)
avgPeakPositionArray.append(avgPeakPosition)
stdPeakPositionArray.append(stdPeakPosition)
numberComponentsArray.append(numberComponents)
sizeLargestComponentArray.append(sizeLargestComponent)
numStudentsArrays[propertyName] = np.array(numStudentsArray)
avgPeakMagnitudeArrays[propertyName] = np.array(avgPeakMagnitudeArray)
stdPeakMagnitudeArrays[propertyName] = np.array(stdPeakMagnitudeArray)
avgTotalInfectedArrays[propertyName] = np.array(avgTotalInfectedArray)
stdTotalInfectedArrays[propertyName] = np.array(stdTotalInfectedArray)
avgPeakPositionArrays[propertyName] = np.array(avgPeakPositionArray)
stdPeakPositionArrays[propertyName] = np.array(stdPeakPositionArray)
numberComponentsArrays[propertyName] = np.array(numberComponentsArray)
sizeLargestComponentArrays[propertyName] = np.array(sizeLargestComponentArray)
studentDistribsByProperty[propertyName] = studentDistribs
# rows = 4+round(len(possibleProperties)/2)
rows = 2
fig, axs = plt.subplots(rows, 3,figsize=(13,3.0*rows))
title = "Simulation for Beta=%.3f, Gamma=%.3f, Pi=%.2f%%, Max class size = %d, in %s"%(beta*24,gamma*24,startProbability*100,maxClassSize,dataSuffix)
if(applyFilterMasks):
title+=" (no CL and no LEC in multi. comp.)"
fig.suptitle(title)
plotIndex = 0;
ax = axs.flat[plotIndex];
ax.set_title("Total number of infected")
for propertyName in possibleProperties:
ax.plot(classesRatios*100,avgTotalInfectedArrays[propertyName],label=propertyName)
high = avgTotalInfectedArrays[propertyName]+stdTotalInfectedArrays[propertyName]
low = avgTotalInfectedArrays[propertyName]-stdTotalInfectedArrays[propertyName]
low[low<0] = 0;
ax.fill_between(classesRatios*100,low,high,alpha=0.15);
ax.set_xlabel("Removed classes (%)");
ax.set_ylabel("Avg. num. of infected");
ax.set_xlim(0,100)
ax.set_ylim(-500,42000)
plotIndex+=1;
ax = axs.flat[plotIndex];
ax.set_title("Infected peak magnitude")
for propertyName in possibleProperties:
ax.plot(classesRatios*100,avgPeakMagnitudeArrays[propertyName],label=propertyName)
high = avgPeakMagnitudeArrays[propertyName]+stdPeakMagnitudeArrays[propertyName]
low = avgPeakMagnitudeArrays[propertyName]-stdPeakMagnitudeArrays[propertyName]
low[low<0] = 0;
ax.fill_between(classesRatios*100,low,high,alpha=0.15);
ax.set_xlabel("Removed classes (%)");
ax.set_ylabel("Avg. num. of infected");
ax.set_xlim(0,100)
ax.set_ylim(-250,20000)
plotIndex+=1;
ax = axs.flat[plotIndex];
ax.set_title("Infected peak position")
for propertyName in possibleProperties:
ax.plot(classesRatios*100,avgPeakPositionArrays[propertyName],label=propertyName)
high = avgPeakPositionArrays[propertyName]+stdPeakPositionArrays[propertyName]
low = avgPeakPositionArrays[propertyName]-stdPeakPositionArrays[propertyName]
low[low<0] = 0;
ax.fill_between(classesRatios*100,low,high,alpha=0.15);
ax.set_xlabel("Removed classes (%)");
ax.set_ylabel("Avg. peak position in days");
ax.set_xlim(0,100)
ax.set_ylim(-2,85)
plotIndex+=1;
ax = axs.flat[plotIndex];
ax.set_title("Total number of infected vs Students")
for propertyName in possibleProperties:
ax.plot(numStudentsArrays[propertyName],avgTotalInfectedArrays[propertyName],label=propertyName)
high = avgTotalInfectedArrays[propertyName]+stdTotalInfectedArrays[propertyName]
low = avgTotalInfectedArrays[propertyName]-stdTotalInfectedArrays[propertyName]
low[low<0] = 0;
ax.fill_between(numStudentsArrays[propertyName],low,high,alpha=0.15);
ax.set_xlabel("Students with at least 1 class");
ax.set_ylabel("Avg. num. of infected");
ax.set_ylim(-500,40000)
ax.set_xlim(0,42000)
ax.legend()
plotIndex+=1;
ax = axs.flat[plotIndex];
ax.set_title("Infected peak magnitude vs Students")
for propertyName in possibleProperties:
ax.plot(numStudentsArrays[propertyName],avgPeakMagnitudeArrays[propertyName],label=propertyName)
high = avgPeakMagnitudeArrays[propertyName]+stdPeakMagnitudeArrays[propertyName]
low = avgPeakMagnitudeArrays[propertyName]-stdPeakMagnitudeArrays[propertyName]
low[low<0] = 0;
ax.fill_between(numStudentsArrays[propertyName],low, high,alpha=0.15);
ax.set_xlabel("Students with at least 1 class");
ax.set_ylabel("Avg. num. of infected");
ax.set_ylim(-250,20000)
ax.set_xlim(0,42000)
plotIndex+=1;
ax = axs.flat[plotIndex];
ax.set_title("Infected peak position vs Students")
for propertyName in possibleProperties:
ax.plot(numStudentsArrays[propertyName],avgPeakPositionArrays[propertyName],label=propertyName)
high = avgPeakPositionArrays[propertyName]+stdPeakPositionArrays[propertyName]
low = avgPeakPositionArrays[propertyName]-stdPeakPositionArrays[propertyName]
low[low<0] = 0;
ax.fill_between(numStudentsArrays[propertyName],avgPeakPositionArrays[propertyName]-stdPeakPositionArrays[propertyName], avgPeakPositionArrays[propertyName]+stdPeakPositionArrays[propertyName],alpha=0.15);
ax.set_xlabel("Students with at least 1 class");
ax.set_ylabel("Avg. peak position in days");
ax.set_ylim(-2,85)
ax.set_xlim(0,42000)
plotIndex+=1;
# ax = axs.flat[plotIndex];
# ax.set_title("Num. of connected components")
# for propertyName in possibleProperties:
# ax.plot(classesRatios*100,numberComponentsArrays[propertyName],label=propertyName)
# ax.set_xlabel("Removed classes (%)");
# ax.set_ylabel("Num. of conn. comp.");
# ax.set_xlim(0,100)
# plotIndex+=1;
# ax = axs.flat[plotIndex];
# ax.set_title("Size of largest component")
# for propertyName in possibleProperties:
# ax.plot(classesRatios*100,sizeLargestComponentArrays[propertyName],label=propertyName)
# ax.set_xlabel("Removed classes (%)");
# ax.set_ylabel("Classes");
# ax.set_xlim(0,100)
# plotIndex+=1;
# ax = axs.flat[plotIndex];
# for propertyName in possibleProperties:
# ax.plot(classesRatios*100,numStudentsArrays[propertyName],label=propertyName)
# ax.set_title("Students with at least 1 class")
# ax.set_xlabel("Removed classes (%)");
# ax.set_ylabel("Num. of students");
# ax.set_xlim(0,100)
# plotIndex+=1;
# for propertyName in possibleProperties:
# ax = axs.flat[plotIndex];
# ax.set_title("Students by level (%s)"%propertyName)
# ax.stackplot(classesRatios*100,list(studentDistribsByProperty[propertyName].values()), labels=list(studentDistribsByProperty[propertyName].keys()))
# ax.set_xlabel("Removed classes (%)");
# ax.set_ylabel("Num. of students");
# ax.set_xlim(0,100)
# if(plotIndex==5):
# ax.legend()
# plotIndex+=1;
plt.tight_layout()
fig.subplots_adjust(top=0.87)
plt.show()
nodeStrength = gclasses.strength();
betweennessWeighted = gclasses.vs["Betweenness Centrality"];
bins = np.logspace(np.log10(1), np.log10(np.max(betweennessWeighted)),23)
plt.figure()
plt.plot(nodeStrength,betweennessWeighted,marker="o",ms=2,linestyle = 'None',alpha=0.050)
plt.xscale("log");
plt.yscale("log");
plt.xlabel("Node Strength");
plt.ylabel("Betweenness Centrality");
plt.show()